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A FORTRAN  PR(X;RAM  TO  PLOT  PHASE  SHIFTS  CP  THE 
TtiREK  DIMENSIONAL  SQUARE  WELL 

ISTROmmON 

Thic  progrut  calculate^;,  end  provides  CalCoatp  plots  of,  phase  shifts  in  redfsat 
or  degrees  for  «echanicsl  scattering  fro*»  a J- 

dimensioftsl  spherically  s^.metrseal  square  veil  potential  of  width  a and  depth  'V, 
for  a da  Sroglte  qua;'.tu»  particle  of  energy  E.  It  exists  in  tvo  \>«tslo«ts,  o*v  of 
which  presides  plots  o'  she  phase  shift  v$  Es® , with  the  (hate  shift  given  In 
units  of  pi  and  in  di^enslionless  unUs«  In  a standard  two  d'MBSional  CalCoatp 
plot.  The  ether  versten  adds  « third  dlMension  to  the  plot,  wit.-  heing  plotted 
v\j  Pa*  and  V.  The  prOijran  penslti  the  first  fis*e  phase  shifts  to  be  calcsilsted 
(l  values  0 thru  V)  and  p.srnits  a detenelnatloo  of  the  mushet  of  e*w«se  bound 
states  the  wel»  will  saps^-trt  vis  Levlnfoo's  ?h«er««. 

1.0  i3E,Vt  IF  ICArii.fi 

1.1  Title 

Caiculattf'-n  and  plc’-tilng  rcutlne  for  deters*! nation  of  pHi-se  shifts  of  the 
spihertcally  sysssetrlcsl  square  well  poterclsl. 

1.2  Entry  points  Both  versions 

PKASHIPT  (2B  and  JO  versiins) 

SPHiR 

MuMarnp*  wSaaeie^  Jwm  tS.  ISTS 


I 


SP«Y!i 


U>'WT  (2D  an<$  3D  vettlocu) 

1.3  Frogrcanins  L«n|;uage 

Unguage:  3^00/3800  rOKOXH 
Soutine  Typ^;  Progr«t,  aubroutiiw 
iperaiing  Syste«;  tStiM  SCOPE  2.1 

1.1  C««si!ter  and  Cont IgMratloc 
CDC  38i>0 

1.5  C»fttiiSjtoT  ov  Pratr&tmer 

J.  L.  feepace,  Code  5216J,  Solid  State  Devices  Srsach,  Elecsrwic*  Tetheoiogy 
DiviSlce^ 

l.h  Contributing  Orgsoiistloo 

KRt  - travel  Research  laboratory,  WasblaftoG,  D.C.  20375 

1.7  Frogratt  Availability 

1.7.1  Subi»ltl;!’.  Progri*  wrtt»nip,  FOtiTliA5i  source  deck,  eaurce  Usttug 

1.8  This  prograae  has  been  used  to  calculate  aad  plot  in  tvo  and  three  diMntiooai 
perspective  the  phase  shifts  ot  the  apherleaHy  ayttaetric  sttuare  well  potential 
ior  square  wells  ts(  depth  0 thru  lOOO  in  danintlottles*  units,  lor  etierfiss  0 

thru  5000  in  dlstensionless  units,  ter  s thru  g-vav*  scatterlnf,.  (V£.*IP1CAI10»> 

1.9  1 June  l'S7.5 

2.0  PV.RF08E 

Althijugh  the  general  propertiea  ef  the  phase  shift  sf  the  three  dliaeasloftal 
sphirically-sfimettrlcal  si^iuare  well  potential*  ar«  destrtbed  in  sN^tst  textbooks 

quancun  siiecKaniea,  spcetfic  graphical  repraseetations  ei  this  function  ate  rare. 

*Shw«''TnT?itri 


2 


this  progrMi  peniits  and  thnf«  diaensiooai  Caicoap  plots  of  s,  p,  d,  f sod  g> 
w«v«  phsfc  «s  a function  of  energy. 

Input  It  hactdled  by  PKASHIFY,  with  the  vull  paruKt^rs  energy  liaits  and 
platting  infcnution  being  read  fro*  data  cards,  the  phaahifts  are  then  calcu* 
lated,  vith  rHASKlFT  calling  SPHVN  and  SPHJN,  routines  which  respectively  ccapute 
.pherical  Seuswnn  and  spherical  Bessel  functiiKt.  Then  TKASKIPT  calls  illSIHAX, 
vhich  rnaputes  the  largest  aod  swllest  pheshifts  calca.lated  if  the  plots  are  to 
be  antoiuticaUy  scaled,  rinally,  the  plotting  routine,  lAYttJT  is  called  and  s 
Cali.ostp  plot  1»  pToducea. 

2.1  Pfobleu  Sackground 

In  aany  field*  physics,  obter^aticsn*  are  described  in  :«nss  of  the  scat- 

tering of  an  electron,  or  other  particle,  on  a target  represented  by  a three 
dlMnsicnal  s^are  well.  I'.xjMaplea  are  found  in  nucles-r  p5'-v'Sic«,  Ir.  eccrsic  phvslcs, 
and  In  solid  State  theory,  Wften,  the  well  appreritt«r**»^  is  only  us\;d 

iniitaUy  to  ohiain  sosse  idea  abfnst  the  ti.ee  of  the  staUere  and  of  *h*  strength 
of  its  potr^ntial.  In  tuny  eases  the  appreKlsMtloo  I*  augiaented  with  % ■ *riety  of 
refinessests.  Keverth.eless,  the  three  dtaensiotial  *;p»are  well  p «?vldt>*  us  vfth 
tlw  prototype  of  the  Ssehavtour  of  a scattering  object.  This  i»  ihs  eah«  In  <{uantu«i 
aecnanit'S,  and  alssu  in  electrcHsagnatic  theory.  Although  tis*  general  properties 
of  the  phase  shift  ef  the  thte*  dlaatncional  spherically'- syasaetrical  square  veU 
potentiai  arc  described  in  woet  textbooks  on  (tvantuae  Mechanics,  stwciflc  g:~fpsi*>cal 
rapeesantations  of  this  fursetloo  are  rat«. 

J.3  CSACK 

3.1  Calli.fsg  S«^f«t»ce  or  Operational  Tresradure 

To  use  gflASRll'^  e*  a prograat,  data  cards  after  she  ’l<ytO  card  are  recited. 

See  section  5,7.  there  are  two  dtffarant  versions,  depending  upon  whether 


piot  dispiays  one  thrau;|h  five  ph£se  sbifca  versus  the  eaergy>vell> 
width  pereJKt.'r  . The  three  di»ensio(ui  plot  dispUjF*  • teiected  phase  shift 
versus  La^  for  up  to  nine  dtffereat  potenclel  depths,  where  the  potential  depths 
are  displayed  on  the  tnird  axis,  which  is  oriented  at  an  anfile  of  about  2CF  frcas 
tht'  abscissa.  ?1<e  well  depths  are  displayed  equally  spaced  aloa|;  the  poteatiai 
axis,  althout$h  they  need  not  be  consecutive. 

1.2  Argujsents,  rsraMeters,  sad  Initial  Cosidivicns 
2-D  Vrrsloft  'od  1*D  Verslor 

{.-swest  energy  value 

r?StL*»K  The  increskent  to  ZDIO  for  each  iteration,  of  which  there  are 
lAST  The  cecal  misber  of  data  points  destred  (B^lLOlii  X - the  KaxiauA 

value  of  tbe  «‘nargyj 
VVI.CAN  the  well  depth 
the  veil  Width 

MAX  The  htglMSt  order  phase  shift  desired  fall  phase  shifts  up  to  MiAX  are 
super itapoeed  on  one  plot! 

'Elii!'  determines  uhether  the  phase  shifts  ere  fiveo  in  decrees  or  radians 
imeWf  supptcsccs  plotting 

ii<&Cfkl£  suppreaset  ordinal  s«!f««oallas.;  if  self-scalinf  is  suppressed,  tss'c- 
fesiA  a.sd  ainleun  values  of  the  ^aa«  i^ift  VS^,  VHIM  he  inserted  prior 
to  statenent  lO'^  in  PiW^HlFT 

KUX  Designated  last  daea  rerd-  the  position  of  KILL  0.1  the  data  tardc  dif- 
fers betwesrt  the  2-0  and  3-0  versions. 

»cr  additional  Infotsution  concerning  these  parsaiietert,  see  Section  3,7. 

3.3  Space  K««iu'red  (Oeclual  Sad  Octal.; 

3.3.1  Unltiue  Storega; 

<■ 


2-D  version: 


PHASHIFT 

12722  octal 

Locations 

MINIMAX 

00230  •• 

SPHJN 

01140  •' 

It 

SPHYN 

00277  •• 

II 

UYOUT 

01332  •• 

•1 

SCALE 

00246  •• 

II 

AXIS 

00536  " 

It 

3-D  version: 

PHASHIFT 

13003  octal 

Locations 

MINIMAX 

00230  " 

II 

SPHJN 

01140  '• 

II 

SPHYN 

00277  •• 

II 

LAYOUT 

n332-  •• 

II 

SCALE 

00246  " 

11 

AXIS 

00536  •• 

II 

3.3.2  Coinnon  Blocks: 

Both  Versions; 

COHMON/IABEL/IL,  1U<5),  ILY(5)  ,ILH(5)  ,ILR(5)  ,1LC(5)  ,ILABEL(5,5) 
CatlON/RANGE/XMAX,XMIN,YMAX,YMlN,IFRST(5)  ,1LST(5)  ,XSIZE 

(In  the  3-D  version  the  subscripts  for  all  variables  and  arrays  in  lABEL  and  RANGE 
are  (9)  and  (9,9)  resp.) 

C0»10N/1/DELTA(1000) .CASTOR(IOOO) .POLLUX (1000) .HERCULES (1000) .AJA 
IX(IOOO) .OLYMPUS (1000) .HADES (1000) .B(IOOO) .DB(IOOO) ,N(1000) ,DN(1000 
2) .BZ(IOOO) .DBZ(IOOO)  .NZ(IOOO) .GAMMA(IOOO) .ENEPGY(IOOO) ,CYCLOPS(100 
30) .CENTAUR(IOOO) .ZEUS (1000) .RHO(IOOO) .DNZ(IOOO) 

CO»10H/A/PHASE(5000),ERGON(5000)  (In  3-D  version  this  is  called  COlMON/3/) 

3.3.3  Temporary  Storage:  None 

3.4  Messages  and  Instructions  to  Operator;  None 

3.5  Error  returns.  Messages  and  Codes 

These  are  all  related  to  the  Bessel  and  Neumann  function  subroutines: 

See  NRL  Computer  Bulletin  29.  section  3.6  for  explanation  of  error  messages  re- 
lating to  overflow  and  underflow  prevention  and  use  of  zero  arg.for  Neumann  func- 
tions. 

3.6  Informative  Messages  to  the  User 

Other  than  the  error  messages  mentioned  in  Section  3.5,  no  informative  mes- 
sages are  printed. 


3.7  Input 

This  section  describes  the  data  required  for  both  the  3-D  and  2-D  versions 
of  Program  Phashlft. 

3.7.1  2-D  Version:  One  card  Is  required  to  control  the  input  of  data.  This  card 

specifies: 


Parameter 

Format 

Columns 

Description 

EZERO 

F10.5 

1-10 

The  smallest  energy  value  desired 
in  the  calculation  of  the  phase- 
shift  vs  energy  plot. 

EPSILON 

F10.5 

11-20 

The  energy  increment  added  to 
EZERO  to  generate  the  next  energy 
data  point. 

VULCAN 

F10.5 

21-30 

The  absolute  value  of  the  well 
depth. 

FERiil 

F10.5 

31-40 

The  well  width. 

11 

45 

The  highest  order  phase  shift 
desired.  For  S,P,D,F,G-waves, 
MAX  = 1,2, 3, 4, 5. 

JELLO 

11 

50 

Controls  units  of  phase  shift. 
If  JELLO  is  1,  phase  shift  is 
expressed  in  degrees,  if  JELLO 
is  blank,  phase  shift  is  calcu- 
lated in  radians. 

LAST 

14 

55-58 

The  total  number  of  data  points 
desired  to  be  calculated,  maxi- 
mum 1000.  The  largest  energy 
value  calculated  is  given  by  tiie 
product  EPSILON  X lAST. 

NOPLOT 

11 

65 

If  a plot  is  not  desired,  enter 
el.  Otherwise  leave  blank. 

NOSCALE 

11 

70 

If  set  = 1,  ordinal  scale  is 
0-1.  If  blank,  autcmatlc  self- 
scaling  obtains. 

KILL 

11 

72 

Zero  except  in  last  date  card  of  set 

3.7.2  3-D  Version: 

Two-plus  cards 

are  required  to  control  the  input  of  data. 

A maximum  of  10  data 

cards  per  set  i 

are  possible. 

The  first  card  specifies: 
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Parameter 

Format 

Columns 

Description 

EZERO 

FIO.7 

1-10 

The  smallest  energy  value  desired 
in  the  calculation  of  the  phase- 
shift  vs  energy  plot. 

EPSILON 

F10.7 

11-10 

The  energy  increment  added  to 
EZERO  to  generate  the  next 
energy  data  point. 

FERHI 

F10.7 

21-30 

The  well  width. 

MAX 

11 

45 

The  highest  order  phase  shift 
desired.  For  S,P,D,F,G-waves, 
MAX  = 1,2, 3, 4, 5. 

JELLO 

11 

50 

Controls  units  of  phase  shift. 
If  set  to  1,  units  are  degrees, 
if  blank,  units  are  radians. 

LAST 

\ 

NERVES 

14 

55-58 

The  total  number  of  data  points 
desired  to  be  calculated,  max- 
mum  1000.  The  largest  energy 
value  calculated  is  given  by 
product  EPSILON  X LAST, 

11 

60 

The  number  of  curves  of  different 
well  depth  desired  to  be  plotted 
(maximum  of  9) . 

NOPLOT 

11 

65 

If  a plot  is  not  desired,  enter 
a 1,  otherwise  leave  blank. 

NOSCALE 

11 

70 

If  set  1,  ordinal  scale  is  0-1. 
If  blank,  automatic  self-scal- 
ing obtains. 

second  and  remaining  cards 

specify; 

VULCAN 

F10.5 

1-10 

The  absolute  value  of  the  well 
depth.  Successive  data  cards 
generally  contain  successively 
larger  well  depths,  up  to  a 
maximum  of  9 total  data  cards 
giving  values  for  VULCAN. 

KILL 

11 

21 

Zero,  except  pn  the  last  data 
card  of  each  successive  data 
set. 

Output 

7 

3.8.1  Printed  output:  &c£h  versions 

The  folleswing  inforsMtian  is  printed  before  the  data  list: 

{!)  PHASE  SmFT  IK  RADIAIB  (The  vord  M<StEES  will  appear  of  SAQJA)€ 

li  that  option  is  selected)  FX»  A qtAimm  HKHANICAL  PARTICLE  SCATTERISG  OF?  A 


’-!)  Sqt’ARii  hELL. 

(2)  Inforaaticti  regarding  the  chosen  Input  parsneters  is  printed  as  follcwa 

'-ILL  !KPTK  = 

WELL  VTOTH  ^ 

LZ£RO  - I The  rsininsa  and  JMxiHua  values  along  the  energy  axis. 

EUST  ■ . 


(3's  tinta  Ltst 
EKTllGl’  si) 
DELTA  i.&y 
miA 
IhlLTA  <2>: 
ffiClIA  i3V. 
aeuA 
K : 

ERGO  : 


7r;o  data  list  contains  IC  variables  as  fcllws; 

JhC'  energy  of  the  lih  data  point 
The  tero  order  phase  shift  (S-»ave) 

The  first  order  phase  shift  (T-vave) 

The  seennd  order  phase  shift  fD»vave) 
the  third  order  phase  shift  (F»wav*e) 
the  fesirth  order  phase  shift  <G-wave) 

Z£fIS(I),  the  protiutt  of  the  sijuare  root  of  the  EKESCVd) 
and  rKSHi^the  well  vldthl 

The  product  S * IliWtl  (this  Is  what  the  actus)  abscisaai 
value  Is  in  the  plots) 


OLYHfifS  1;  The  nuaeratox  of  tho  Arctan  fuoctioo  for  1 ) data  point 

<printed  for  dlagncstic  pjrpc-ses  In  checking  the  accuracy 
of  the  phase  shift  caloulaiion)  tee  sect.  4,0. 


HABIS  1 • the  dsnoakinator  of  the  Afctaa  function  described  above. 


3,8.2  Flotted  O.Jtput 

1-n  version;  Sasiples  of  the  plotted  output  ate  given  in  fii^res  A,  S,  and 
fe.  Each  of  these  fi^^at**  shows  the  S,  P,  9,  F,  and  G*wave  ;^te  shifts  supesrlw- 
;mt6  on  I'iw  satae  plot  for  cowpatative  purposes.  The  ordinal  value*  glvi*  the 
jshsse  shift*  in  unit*  of  pi,  while  the  abscitsal  values  are  product*  of  energy 


8 


X well  Width‘S.  The  cur  have  b«en  «afc«Lllsh£d  for  Illustrative  {Kxrposea, 

the  actual  curves  are  a.,  displayed  at  solid  lines,  and  luat  be  identified  indi- 
vidually either  by  plotting;  consecutive  displays  uorking  up  froa  S to  C-wavea, 
or  b«  recourse  to  the  printed  output.  The  interpretation  of  the  broken  curves 
displayed  In  figure  b is  discussed  in  section  S.O  (KEHAJiKS);  it  is  a result  of 
the  i,rvli8tlnguishabllity  of  the  sultipies  of  out  to  the  coaiputer. 

3-D  version;  A saaele  o'  the  plotted  output  is  given  in  f&gura  3.  xt  is 
eafee  111  shed  {Stly  by  the  nunerals  identifylni;  V-values  necr  each  curve,  these  are 
not  displayed  on  the  CalCoap  plot,  otherwise  the  figure  la  essentially  uoretouched. 

It  should  be  noted  that  it  is  not  necessary  to  have  sequential  V-valuet;  it  should 
also  he  twted  that  the  eaxlsasu  point  of  each  individual  V curve  is  displayed  be» 
tween  0 and  I in  units  of  pi,  atvl  in  no  case  Is  greater  than  1.  la  casca  whare 
the  curw  would  have  penetcatee  tl««  region  between  I and  2,  the  coatputtr  dlsplaya 
the  peak  of  the  cur\‘e  reduced  to  the  Interval  0 to  1,  as  in  curve  V 20. 

3.9  Fonastt 

In  the  2-D  veraxon,  Che  data  are  printed  in  FlO.3  t'oruat;  in  the  3-3  va^’sioa, 
the  data  ate  printed  in  ?U.7  forwat. 

3.10  Fvtetnal  Routines  and  Sya&N^ls 

both  vetsi.v';si»-.  ftikSSm  tails  SHIFK 

&ms 
umimi 
LAVOlt 
STCPVUff  J 

UUK2  ] 

Si)g,tv  ) -SiSTEH  UB&iS¥ 

3. 1 1 TlMng 

2-D  version'  To  produce  fituraa  such  aa  fii^ces  4,  3 and  b,  i.a.,  foe  lOGO 
data  points  and  calculations  of  S thru  C wave  pHaae  shifts,  the  procsstingi  tins 
U generally  less  than  1.3  ainutes.  The  plotting  tiua  generally  runs  about  9 «Uutaa. 

§ 


( _ 


3-D  Version:  To  produce  a figure  such  as  figure  3 generally  takes  aboi\t  4 

minutes  of  computer  time  plus  about  4 minutes  of  CalComp  time. 

3.12  Accuracy 

The  accuracy  of  the  calculations  was  verified  using  an  electronic  calculator 
for  several  points  on  the  S-wave  curves.  The  accuracy  of  the  plots  was  verified 
by  comparison  with  curves  for  S,  F,  and  D-wave  phase  shifts  presented  for  V = 

57.5  over  an  Ea^  Interval  of  0 to  10,  which  were  presented  In  Mott  and  Massey^ . 

3.13  Cautions  to  Users 

Experience  has  shown  that  the  main  sources  of  error  In  using  these  routines 
lies  In  the  confusion  which  may  arise  In  identifying  tne  individual  phase  shift 
curves  on  the  Calcomp  plot,  especially  at  the  hl^  energy  end  of  the  curves. 

A second  source  of  confusion  may  arise  from  the  display  of  broken  curves  In 
figures  3 and  6 due  to  the  compression  by  the  computer  of  multiples  of  pi  Into  the 
first  "zone"  of  zero  to  pi.  This  reduction  is  demonstrated  In  figures  2a  and  2b, 
the  uiureduced  and  reduced  phase  shifts  as  applied  to  figure  5a,  page  33,  in  Mott 
and  Massey.  A method  of  determining  the  ordinal  Intercept, so  as  to  be  able  to 
apply  Levinson's  theorem^ which  relates  the  ordinal  Intercept  to  the  nund>er  of 
bound  states  the  well  may  support,  Is  discussed  In  section  8.o,  REMARKS. 

3.14  Program  Deck  Structure 

2-D  version:  JJOB  card 

JEQDIP,  13  = PL 
JFTM  card 

Program  Subroutines:  PHASHIFT 

SPHJM 

SPHYM 

MINIMAX 

lAYOUT  (with  COMMON  blocks  and  DIMENSION  statements 
SCOPE  subscripted  5) 

JLOAD  card 
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^RUM  card 


data  cards  (see  sect.  3,7.1)  each  data  card  will  generate  a series  of  sup- 
erimposed phase  shift  curves  for  a given  well  depth;  each  successive  data  card 
will  generate  another  plot. 

E0F 

3-D  version:  ^JOB  card 

2 DEMAND,  53776b 
JEQUIP,  13  = PL 
JFTN 

Program  subroutines:  PHASHIFT 

SPHJN 

SPHYN 

MINIMAX 

LAYOUT  (with  COMMON  blocks  and  DIMENSION  statements 
SCOPE  subscripted  9) 

’BANK,  (0«,/3/) 

’LOAD  card 

’RDN  card 

data  cards  (see  section  3.7.2) 

E0F 
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4.0  METHOD  OR  ALGORITHM 

I.  A detailed  description  of  the  significance  of  certain  features  of  the 
phase  shift  carves,  as  well  as  a guide  to  the  literature  Is  available  In  reference 
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(2).  Th«  MtheiMtlcal  expreitlOQ  for  ch£  pbtu  •hifc*  «nd  siapllfyia^  ass<Biptioa« 
resardio^  the  cotutenta  are  discuaseo  belcr*, 

4.1  DeacripCloo  of  the  Phaae  Shift  &^clooe 

Ve  ROW  consider  the  ideas  of  the  partial  wove  Method.  Stsppoae  we  have  a 
solution  inside  a spherical  well  of  radius  a and  depth  9.  Sec  figure  1.  this 
solution  has  to  ststch  the  wavefuacticoa  outside  the  well  at  r » ■,  At  this  point 
the  waveftmetion  should  have  the  sa»e  value  (coatinuitp  of  the  wavefunctioo)  and 
the  sane  derivative  (coatiouity  of  the  derivative  of  the  wave  fttoction}.  The 
continuity  of  the  derivative  of  the  wavcfunctioa  is  retired  because  the  nosbentua 
of  the  particle  should  be  defined  at  every  point  in  spaco.  There  are  two  store 
boundary  conditions:  1;  the  wavefunction  inside  the  sphere  should  be  tero  at  the 
ortiiin,  at^  2)  the  wavefunction  outside  the  sphere  should  Save  s fined  aaspUtude. 
The  value  of  this  saplitude  is  of  no  ie^rtaoce  when  we  csUulate  the  p^se  shift. 

first  boundary  condition  excludes  heunaan  functions  inside  the  sphere.  As  a 
result  of  the  spherical  synssetsy  of  the  potential  well,  the  angular  depenlcnce 
‘.Rside  iind  outside  should  be  identical:  the  i*valu«  outside  the  sphers. 

We  intrewkice  the  following  notation,  k ^vE.  W use  units  si  h ^ 1 alsd 
■t  V where  V ' the  well  depth,  the  sbsotuce  value  of  the  negative  potential 

Inside  the  sphere.  Ue  use  the  staadlard  'location  for  the  spherical  lasset  and 
Hewasisn  functions  (coapare  idi3  handbook The  radius  of  the  sphere  is  r - s. 
With  this  nocation  we  have,  far  the  radial  part  of  the  wevefunctina  Sasida  the 
epftere : 

r<a:  *<e>  * Aj|f>r) 

and  for  the  wavafunctlon  cueside  the  spbers: 

r?«:  s(r)  «•  B}|fkr)  a Oa^(kr),  * 1 

The  convention  Is  nomelly  to  replsca  the  two  conatante  by  the  sine  and  cosine 
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of  an  an^ie  at  follows . 

Having  ettablisKe<3  there  two  wavefuactlona,  the  aatchlng  coauHcico  givaa  the 
following  relation: 

i^(ia)  - c&a  - Bin  5^n^(lui) 

- k (coa  $^j^(ka)  - alt*  5^0|(ka)) 

We  take  the  ratio  between  these  two  e<uatiotia.  oaltcing  the  subscript 

j‘(i>a)  kj’fSial  •>  kun  5,  p'(ka) 
j(>a)  J(ka)  - tan  { nOta) 

and  solve  for  tan  the  result  can  be  vrlttcn  as 
* tan-1 

^ka  n'cks)  - n(ka)- 

where  y is  given  by 

V ja  }‘<*a)/J(aa) 

These  are  tire  e^atiohs  of  interest  which  ere  used  to  ccaiiptita  the  phase  shift 
b.O  Lktting  of  pro^raa)  ~ tte  pagN  2-5  «<  «n4 
b.O  OOWARlSOMi  « see  section  5,U 
7.0  rest  HEtKOO  ASD 

S-D  version:  To  g^enerate  fig^sree  it,  S and  b during  oua  eoahniter  run,  set 


1st  date  card 

itki  data  card 

3rd  da 

,0<KKU 

.GOOsJl 

.00001 

ETSllft^  0.5 

e.s 

0.5 

W‘5J[aii}i  : 1.0 

10. 0 

100.0 

IT.«ll  1.0 

1.0 

l.O 

5 

5 

JEIJU)  - b 

b 

b 

usr  IOC 

IGO 

100 

lioruor  : b 

b 

b 

ItSOALE  - 1 

1 

1 

klU  0 

0 

1 

the  partial  printed  output  item  the  second  data  card  follows;  (See  Figure  A) 
3'D  version;  to  generate  figure  I,  set 
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1st  data  card 


E2EE0 

EPSILON 

FKSMl 

HUC 

?r!xo 

LAST 

Ka'K%’£S 

ViTtSCAJ  E 

laid  4ats  .atd 

’ ” 

*■  th  ^ 

9 „ \ 

ft*' 

^ih  ■*  " 

i Jlh  " 


.OOOOOOl 

.014 

1. 

5 

b 

1000 

9 

b 

1 

VULCAN 


1 

2 

3 

4 
10 
18 

19 

20 
21 


K.1LL 


b 

b 

h 

b 

b 

b 

b 

b 

I 


\ viri'.il  listing  of  prtaat«d  output  fer  figure  3 follows:  (s««  Figure  B) 

?i,  t»i  t^ai  the  pisasr  shifts  arc  sUfftnad  only  to  within  an  integral  auUiple 

:!  pi.  U Is  CTsStoaaty  to  Ua  *«  > 0.  Since  the  ewsputer  delivers  values 

« -v-^  * 

*4  phasi*  shift  only  between  0 end  1 in  units  of  pi  radians,  the  t^ueaeioo  arises 
s»  how  to  interpret  these  plots  in  li|?vt  of  leelnsoe's  thcctea,  Aerording  to 
l^  .tns^vn's  thecreie.  she  value  ot  the  phase  shift  at  tero  en«r8,y.  in  units  of  pi 
tadiant,  tieVis  the  nundHtr  of  |-«av«  bowwJ  state#  that  the  well  can  »upt>ort,  ex* 
«apt  at  certain  transitional  strengtha*.  If  Multiples  of  pi  ar»  auppressed  in  the 
plf?t,  hnr!  ran  c«»«  detemine  the  taatber  of  Ivound  itates?  It  will  be  noted  in  figure 
b for  the  deep  well,  that  there  are  a nuaber  of  diiconitnuities  in  phae*.  this 
is  a result  of  she  suppressioo  of  tatltiples  of  n pi.  If  cne  plots  figure  B out 
so  energies  greater  than  400,  it  will  be  fauad  that  only  one  aura  sat  of  dlscon* 
tituitiea  occurs  before  the  phase  shift  curves  assywptotically  approach  sere.  If 
one  assifo#  a negative  value  to  Juags  up,  aad  a poaitive  value  to  JjuMpa  down,  the 
aigebraic  nunhar  of  phasa  juags  will  be  agjual  to  the  nuabar  of  l*w>«v«  bautwi  atates 
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9.0  4CKMVL£Z>S£3fiEirr 

The  suthor  ic  Inicbted  to  Pcofeator  P.  M.  E.  Heijer,  who  •ofseeted  thle  wock 
•ad  to  Hr.  S.  Bteaacr,  who  coacrlbuted  •ubroutiaet  LAYOiTT,  SCALE  tail  ASIS,  to  Kre 
J.  KcaoQ,  who  cootributed  eubroutln««  STHIK  SPSTK,  end  to  Dr.  A.  H.  Lowr«;r. 
for  helpful  discucsioo*. 
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2-D  Version 





PROGRAM  PHASHIFT 
C Two  DIMENSIONAL  VERSION 

c This  program  calculates  phase  shifts  in  radians  or  degrees  vs  e*a**2 
C OR  K OR  E IN  DIMFNSIONLESS  UNITS  FOR  ThE  OUANTUM  MECHANICAL  SCATTERING 
C OF  A PARTICLE  OF  FNFRGY  E FROM  A 3-D  SqUaRE  ^ELL  OF  WIDTh  A AND  DEPTH 
f -V,  The  input  parameters  are  « ii.  e^ero  (The  lowest  energy  value)  co 

COLS  l-lO.I?).  EPSILONIThE  energy  INCREMENT)  COlS11-20,  (3).  /ULCAN 
r TmF  absolute  VALUE  OF  ThF  WELL  DEPTH)  COLS  21-30 

r (a),  eermi  (The  well  width)  cols  3i-ao,  (5)»  max  (The 

C HIGHEST  order  PHASF  SHIFT  DESIRED.  FOR  S.P,D*F)  MAX<1«2<3>A.  COL  A5 
C (ft).  X!LL(N0N2ER0  ONLY  ON  LAST  DATA  CARD)  COL  72.  (7).  LASTIThE  TOT 

C At  number  of  POINTS  DESIRED  TO  BE  CALCULATED).  MAX.  lOOO.  COLS  55-58 
C (p).  JELLO  (IF  PRINTOUT  OF  PHASE  SHIFTS  IN  DEGREES  IS  DESIRED*  ENTER 
C 1 IN  COL  50.  IF  RADIANS  DESIRED.  LEAVE- COL  50  BLANK)  COL50 
C IF  NO  plot  is  desired.  ENTER  A 1 IN  COL.  65.  IF  ORDINAL  SELF-SCALING 
C DESIRED,  enter  A 1 IN  COL  70.  AND  INSERT  VALUES  OF  YHIN  AND  YMAX  BEFO 
C T CALL. 

C0MM0N/LABFL/IL.ILX(5) »IlY(5).ILH(5) .ILR(5).ILC(5) »ILABLE(5.5) 
COMMON  /range/  XMAX,XMIf!*YMAX*YMlN.IFRST (5) *IlST(5)  .XSIKE 

real  ilx.ily.ilh.ilr 
real  n 

REAL  NZ 

COMMON/1/  DELTA! lOOU) . CASTOR! lOoO) .POLLUX ( 1000 ) .HERCULES ( 1 000 ) .AJA 
1X( ICOO) .OLYMPUS (1000) .HADES! IOCl  ) .B( 1000  » .DB( 1000) *N( 1000) .DN( 1000 
?) .B2( 1 J('0) .OPZ! lOO'M  .NZ(lCOO) .GAMMA! 1000 ) .ENERGY  I lOCO ) .CYCLOPS! lOO 
PC)  .CENTAUR!  ICC '))  .ZEUS ( ICOC ) »RHO!  IC-DO ) .DNZ!  lOCO) 
dimension  DELTaI  ! 10>-C  ) .DFLTA2(  lOoO)  »DELTA3!  1000)  .DELTA 

1 A( IOOl ) .W( 5) 

DIMENSION  EPGO(lOOO) 

common/ A/  PHASE (5000) .ERGON (5000) 

0 PFAD  (60,10)  EZER0«FPSIL0N.VULCAN*FERMI .MAX.JEUO.LAST.NOPLOT .NOSC 

ialf.kill 

1 0 format (AIFIO, 5). ax. II, AX»I1.aX.I<,.6X. Il.AX.Il.lX.il) 

IFPST! 1 ).l 
IERST(2).1-.  last 
IFRST(3)»1+  7*LAST 
IFRSTIa)*!-.  3*LAST 
1FRST!5)  = 1-.  A*LAST 
ilST(i).  last 

ILST(2).2*LA5T 
ILST (3)«3*LAST 
ILST ! A)*A*LAST 
ILST!5)»5*LA5T 
IE  (JELLO. EO.l)  GO  TO  2 

1 print  6 

GO  TO  3 

? print  7 

3 print  17, VULCAN. FERMI .FZERO.XMAX 

PRINT  8 

6 FORMAT  !1X,*PHASE  SHIFTS  IN  RADIANS  FOR  A QUANTUM  MECHANICAL  PARTI 
KLE  SCATTERING  OFF  A 3-D  SQUARE  WELL.*/) 

7 format  !iX,»PHASE  SHIFTS  IN  DEGREES  FOR  A QUANTUM  MECHANICAL  PARTI 
KLE  scattering  OFF  A 3-D  SQUARE  WELL.*/) 

A FCRMAT(5X.*ENERGY! I ) *<,X«DELTa ‘ 0) *aX*DELT A ! 1 ) *aX*DELTa ! ? ) *4X*DELTA ! 
^3)«(,X»DELTa(4)*(.X»  r.  *4X*  ergo  *<,X«0LYMPUS1*4X*HADES1*/) 

!•>  FORMAT! IX. *WELL  PEPTh«  •FIO.5,*  WELL  WidTh*  *F10.5.*  EZeRO»  * 

1 F10.5.*  ELAST.  *F1C.5/) 

ICsMAX-fl 

'■NERGY(1).FZER0 


100 

200 

300 

400 

500 

600 

700 

BOO 

900 

1000 

1100 

1200 

1300 

1400 

1500 

1600 

1700 

1800 

1900 

2000 

2100 

2200 

2300 

2400 

2500 

2600 

2700 

2800 

2900 

3000 

31C0 

3200 

3300 

3400 

3500 

3600 

37P0 

38C0 

3900 

4000 

4100 

4200 

4300 

4400 

4500 

4600 

4700 

4800 

4900 

5000 

5100 

5200 

5300 

54OO 

55OO 

5600 

5700 

5800 

5900 
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no  >)l  t-l.LAST 

FROOI  n«eNFROY(  n»(FERMl  »•*? 

U FNERGYI  I+I  )»FPSILON'*-FNFRGYI  I » 
lA  fYCUOPS( n.SORT(ENFRGY( n I 

?n  CFNTaUR(  1 )«SORT  (FNFRfiYt  I )*VULrAN) 

?A  ^FUS(  I )«CYrLOPS(  n*FERMI 
^EUSS»^FUS(  I ) 

77  IF(2f0SS,F0,01  78*30 
■*0  RHOl  n»CENTAUR(  n»FFRMI 
RHOO»RHO( I ) 

7*  CALL  SPHJNl  PHOOtIC  «B) 

CALL  GPHYNlRHOOtK.Nj 
DB>  n=-ni?) 

OM  n»-N(  2) 
no  3A  J!«?.MAX 

np<jn.(ji*p(ji-i)-<ji  + i )*Rt  ji+i n/i2*Ji*i ) 
nNi  jn  = i ji»Mt  ji-n-(  ji4i  )»N(  ji+i)  i/(2*ji+i) 

3A  continue 

CALL  SPHJN( ZFUSS.K ipZ) 

CALL  SPHYNtZFUSStK*NZ) 

PrZ( 1 )=-oZ(2) 
nN2H)r-NZ(2» 

77  no  3F  KI»?.MAX 

ORZUI  l = tKI*PZlKI-l  l-Ul  + l »«RZ(|Cl4l  ) J/l2»Kt4l  ) 

ONZIicI  1 X (KI»nZ  (If  I-l  + i )*NZ(Kl4in/I  2*IC  1 + 1 ) 

3fl  continue 

AO  no  7ft  Llxl .MAX 

Aa  r.A"MAlLI  )xCFMTAUS(  n *np(LIl/B(LI) 

AA  CASTOPd  I 1«CYCL0PS(  n*OPZ(LH 
ftO  POLLUX(Ll  ;xC.AMMA(Ln*PZ(Ll  > 

At  mfpculft.(li  ) xcyclopsi  I »»0NZ(Ln 

ftO  AJAXII  n=rAMMA(LIt*N2lLI  1 

ft.;  Olympus iL I ixCastorili  )-pollux'.li  i 

7P  HADFStLl )xHrPCULFS(Ll I-AJAx(lI  1 
7t  U(LI1  xaTaN2<0LYMPUSIL1) »HAnFSlLn ) 

IF(W(lI t .LT.Ul  U(L1 )=W(U1 U3.14159 
IF(  JFLLO.FO.C'J  go  to  735 
Wd.l )=W(LI )»57, 79578 
IF  ( JFLLO.EO.I  IfiO  TO  76 

7Ar  W ( L I )=W(  L I ) /+,  1 Al«iO 

7ft  CONTINUF 
77  GO  TO  80 
7fl  no  79  LI«1»MAX 
UlLllxO. 

7q  CONTINUE 

OLYMPUS! 1 ) tO, 

HADES! 1 )xO, 

8 0 print  R a,  ENERGY!  I ) .W  ( 1 ) .W  ( ? ) .W  ! 7 ) »W  ( A ) .V.  I 5 ) .ZEUS!  I ) .ERGO!  I i .OL 
IYmpuS! n .hAOES!  1 ) 

8t  format  t lO  1 7X.f10,5 ) ) 

DFLTaI  1 ) tV»  ! 1 ) 

DFLTaI ( I) xWl 2) 
delta?! I )=W(7) 
nFLTA3( I )=W(A) 

OFLTAA ( I 1 tW( 5) 

PJ  CONTINUE 
8ft  print  87 
87  format  (///) 

IFInOPLOT.FQ.1 ) GO  To  95 


6000 

6100 

6200 

6300 

6AOO 

6500 

6600 

6700 

6800 

6900 

7000 

7100 

7200 

7300 

7A00 

7500 

7600 

7700 

7800 

7900 

8000 

8100 

8200 

8300 

8A00 

8500 

8600 

8700 

8800 

8900 

9000 

9100 

9200 

9300 

9A00 

9500 

9600 

9700 

9800 

9900 

10000 

10100 

10200 

10300 

lOAOO 

10500 

10600 

10700 

10800 

10900 

11000 

11100 

11200 

11300 

llAOO 

11500 

11600 

11700 

11800 


m 


IF  IN0SCALE.E0.1IG0  To  99  11900 

CALL  MINIKAX(OELTA.VMIN.YMAX.LAST)  12000 

99  XSIZE«lO,  12100 

XMIN*ERG0«1I  12200 

XMAX«ERGO»LASTi  12300 

IF  SElF-SCALING  not  USE09  THEN  YMIN  AND  YMAX  VALUES  MUST  BE  SUPPLIED.  12*00 

YMIN«0.  12500 

YMAX>1.  12600 

IF  DESIRED.  OIFFFRFNT  ORDER  PHASE  SHIFTS  MAY  BE  PLOTTED  ON  THE  SAME  12700 


graph.  use  NOSCALF  To  bypass.  PUT  A 1 

IN  COL  70.  supply  YMIN.YMAX. 

12800 

NLAST=S*lAST 

12900 

j.last 

13000 

DO  90  M»1.J 

13100 

phaseim).oelTa(m» 

13200 

Phase  <;^+ji»deltai  iMj 

13300 

PHASF«M+i»*JI*0ELTA2«M) 

13*00 

Phase <m+3*j>*dflt A3 (Mj 

13500 

PHASE(VI-.**J>»OrLTA*IMI 

13600 

ERG0NIM)sFRGO(MI 

13700 

FRGONlM+jj.ffiCiOlMl 

1 3000 

FRG0N<M*3»J1«ERG0(M) 

13900 

FRG0N(M+3*JI«FRG0(M1 

1*000 

FRGON(M+A«jj=EPGOtMj 

1*100 

90  continue 

1*200 

CALL  LAYOUT (ergon .Phase. NLAST.6HE*A**2»7H-DELTA-. 6.7. 51 

1*300 

PRINT  luo 

1**00 

luo  FORMAT!///) 

1*500 

9a  IFUILL.FO.O)  go  TO  9 

1*600 

CALL  STOPPLCT 

1*700 

199  END 

14800 

subroutine  fMNIMAXI APRAY.YMIN.YmaX. 

LAST) 

14900 

dimension  APPaY(IOOO) 

U-000 

This  subroutinf  finds  poth  the  largfst 

AND  smallest  numbers  IN  ARRAY. 

15100 

largest  number 

15200 

J=1 

15300 

15*00 

in  DO  30  KrjJ.LAST 

15500 

IE  ( ARRAY  u ) ,GF.  ARRAY  (KM  30.20 

15600 

20  J = K 

15700 

if(<,eo.last)  go  to  *u 

15800 

GO  TO  SO 

15900 

■>0  CONTINUE 

16000 

40  YMAX=ARRAY(J) 

16100 

print  4S.JJ 

16200 

4K  format  (1X«*JJ=*I3> 

16300 

GO  TO  6»’ 

16*00 

SO  JJsJ.fl 

16500 

GO  TO  10 

16600 

S’.aILFST  number 

16700 

60  Lsl 

16B0U 

LL=lfL 

16900 

lOo  DO  300  M.LL.LAST 

17000 

IE(  ARRAYILI  .LF.ARRAY(Mn300.?00 

17100 

20o  L=M 

17200 

TF(M.F0,LAST)G0  to  400 

17300 

GO  TO  SOO 

17*00 

300  continue 

17500 

400  YmINsARRAVILI 

17600 

print  4S0.LL 

17700 

25 


<»S0  format  (1X*«LL»«I3) 

00  TO  600 
500  LL-L+1 

GO  To  lOO 
60o  return 
END 

SURROUTINF  SPhJN t X «N t APR ) 

NOL  BESN  BESSEL  FUNCTIONS  OF  THE  FIRST  AND  SECOND  KINDS 

dimension  ARRtn 

RX»1,0E+307 

NF»1 

NN»N+1 

IF(X.EO.O.O)  GO  To  A 
IFIX.LE.O,5)GO  To  6 

C FOR  ORDERS  LESS  THAN  ARGUMENT.  AND  ARGUMENT  GREATER  THAN  0.5 

ARRI l)»SlNlXj/X 
ARR(?).ARR(1>/X-C0S(X)/X 
NM«X-2 

IF<N.LT,X)MM«N-1 

IF<NM.LT,0)N'.1»0 
DO  16  I«1 »NM 
FN=I+I 

ARR( !♦? 1 r ( FNfl .0) *ARR( I+l »/X-ARR( I ) 

C CHECK  FOR  UNDERFLOW 

IF(  ARRI  I+?)  .NE.O.OGO  to  16 
ARRI I+l )=ARR( I+lt*RX 

ARRI  I-.2)  = IFN-.1,0)*ARRI  1 + 1 l/X-ARRI  1)»RX 

print  1A.  X.I.N 

16  CONTINUF 

IFIN.LT.XIRFTURN 

H0LD=ARRINm+2) 

r FOR  ORDERS  GREATER  TriAM  OR  LOUAL  To  ARGUMENT.  AND  ARGUMENT 

C GREATER  THAN  0.5 

1«  ALPHA=1.0 
Mo-2-MM 
NOTe=L 
1=N+a7 

IFIX/N.GT.  >,9) I=l+X/lo.O 

ARRI 1+2) =ALPHA 

ARRI  1+I)s0,0 

I=-I-l 

no  2 KK«I.M 

k=-kk 

M00=K+1 

TkP1=K+K+1 

ARRIK)=T<Pi*ARRIC+1)/X-aRRIK+2) 

IF  EXPONENT  fault  1P»2 

1 a TfMpY  = ARRIK )*ARRIK  + 1 ) /RX*TkP1 /X-ARRIK+2)/RX 
ARRIKI=ARRIK)*RX 

if  exponent  fault  2A.2 

2a  ARRIK)=TfMPY 

ARRlK+1 )=ARRIK+1 ) /RX 
PRINT  n.X.N.NOO 
2 CONTINUF 
1?  C=H0LD/ARRIK1 

ifIc.ne.o.oigo  To  ia 

H0LD*H0LD*RX 
N0Tf=NOTF+1 
GO  To  12 


17000 

17900 

18000 

18100 

18200 

18300 
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19700 
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20000 

20100 

20200 

20300 

20A00 

2050C 

20600 

20700 

20800 

20900 

21000 

21100 

21200 

21300 

21*00 

21500 

21600 

21700 

21800 

21900 

22000 

22100 

22200 

22300 

22*00 

22500 

22600 

22700 

22800 

22900 

23000 

23100 

23200 

23300 

23*00 

23500 

23600 


U IFlNOTE.FO.OGO  To  17 
MULT»107*NOTE 
NUL»L-1 

PRINT  ?1  iXtN*MULtMULT 

?1  format  (-<X»F0R  ARGUMENT  OF  *F17.lO*  ANO  ORDERS  OF  *IA*  To  *I4* 
1 *SPH.im  exponent  was  INCRFASED  RY  *16) 

17  on  s I»KiNN 

ARR(  I )»A7R(  n*C 
CONTINUE 

return 

/,  ARR(1)b1,0 
DO  C'  L = ?iNN 
ARR(l)=U,0 

S cont:;mue 
return 

c series  for  argument  less  Than  or  foual  to  0,5 

FP5*S,0F-n 
CH=cs=nN«r>ivai.o 
no  lO  L*1»NN 

ITeXP=NFXPA»MEXPA=LFXPA=0 

NO=L-l 

IF(L.EO.l)  GO  TO  7 
CS»CS»X/DN 

C CHECK  FOR  underflow 

IFIcS.NF.O.U)  go  to  7 

C SCALE  FOR  THIS  N AND  ALL  FOLLOWING  N 

CS=CH*RX*X/ON 
NFXPA*NFXPA+1 
7 DN=nN+?.0 
CH  = CS 
XKK=0.0 

A=CH2=F0PTT=1.0/0IV 
Tn=2.0*NO 
FN=-0.5*X»X 
A XKK=XKK+1.0 

F0=XKK» ( TM+XKK+XKK+1 .0 
A=A*FN/FD 

C CHECK  FOR  UNDERFLOW 

IF< A.NE.T.UIGO  to  25 

c Scale  for  this  n and 

ArCH2*RX*FM/FD 
DIV=OIV/rX 
MEXPA=MEXPA+1 
F0RTT=F0RTT*RX 
25  CH2»A 

TT=A+F0RTT 

IF(TT.E0.0,0,0R.ABS(A/TT).LT,EPS)  go  to  9 
FORTTrTT 
GO  TO  8 

o ARR(L)<=CS*TT 
C CHECK  FOR  UMDFRFLOW 

IFI ARRIL) ,NF.O.O)GO  TO  51 

c Scale  for  this  n and  all  following  n 

ARH(L)=CS*RX*TT 
DIV*DIV/RX 
LEXPA*LFXPA+1 

51  ITEXP=NEXPA+MEXPA-cLFXPA 
IF(  ITEXP.FO.w)GO  to  10 
IF(  ITEXP.KO.DGO  TO  29 
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23800 
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print  NFXPA»MFXpA.LFXPfl 

■>0  FORMATMX.*  f.CALE  FRROR  •Il’XtMOn 
?Q  print  n*  XtNO*N 

l->  format  ( iX««FOR  AROUHENT  of  ••F17.10**  and  ORDFRb  PF  **lA»*  To  ** 

\ IAt»  fiPHJN  FXPONFNT  WAT.  INCPFASFO  PY  107*1 

1C  continue 

FND 

SUBROUTjnF  SPhYnI X *N«ARR) 

niMFNSlOM  ARRI 1 ) 

RX«1  ,0F*1U7 
IFlX.FO.'^.I'^O 
t,  print  a 

t format  ( IX «*nn  MOT  ur.F  zfpo  apoumcnt  for  sphyn*i 
return 

r find  YO  aNO  Y1 

■»  ARR  ( 1 ) =-fOr>(  X) /X 

ARR(?)tAPR(1)/X-5U4(X)/X 
Nm1=N-1 
CO  ? K.  = l tN»'l 
TxPl=Y+X+l 

ARR(K  + ?)=TifPl*ARR(F4.11/X-APRIK) 
r CHFCX  FOU  nvEPFl.nW 

IF  EXPONENT  fault  fl»2 

0 APR(X+1 )=ARR(K*1 ) /RX 

ARRIF  + .?)  = TfP1«  APRIF  + 1 t/X-ARF  IX)/PX 
print  7«X«»’«N 

1 format  ( IX  **F('R  APr.UtEfiT  OF  »tF17.lO.»  AND  ORDFRR  OF  **!♦**  TO  ** 
\ IA»»  0PhY»!  FXPON^^'.T  WAS  ^FCPFASEO  RY  1u7*  I 

7 CONTINUE 


FMO 

SUPOOUTINF  LAYOUT  ( XP FFF  • YP.FF R *NN ♦ X(>C0 * YoCp .I^NX  *NNY. NUMBER  ) 
COMmon/LAPFL/ IL* ILX ( S I • ILT (S) • ILH( A I •ILRIf  ) VlLC  ) * ILARLE  t A*S) 
common  /RAnOF/  XmaX«XA;IN  .YMAX.YM)N»IFRS>»S>  ♦W.‘o*T-t5  ) *XSUE 
real  ILX.ILY.ILHtILR  . , 

OtVFNSlOM  XPCCISl  •YPCO(A), 

dimension  XPFEMNMI  .YPFFRlnf;!  .PLTAR^RAYT25a) 
oaTa  II STRT=1 1 • * j * 

initialization^  ^ 

N*NN 


NX=NNX  . . « - 

NYsNNY  • 

IF  I XSIZf.pO.U.UI  X5VTtF=lO. 

ysize'Iu. 
nsvf=n  y' 


n«“»psf  iaij 

IF  (NUMRFR.FO.O)  NUMPfP-1 
IF  I tSTRT.FO.2 1 00  TO  ? 

CALL  PL0TS(PLTARRAY,2AAt:i) 
iStaTi? 


continue 


lapel  X AXIS 


'‘X*  (X-AX-XfIN)/XSIZE 
IF  OX.SF.'  1 00  TO  1 

CALL  SCAlF  I XBFFRtN iXSIZF *XMINtOXtl ) 
'V'  TO  7 
1 COM1»cUF 
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o r\ 


»riiniwn  IT  >Hii 


r print  «>'iS'i*XBFFR(  II  •YRFfRfl  ) »1 

format  (?F2t).lO.IlU) 

XRFFR(n=  (XRFFRm-XMINI/PX 
7 fONTlNUF 


FAU.  f>CALF  (XPFFR*N»XSl2E*XMIN«OXtl  I 

/FRO- XSI 2E-2ERO 
tF  (NSVF) 

* XMtN  = X3l2F*r)X-ABSFtXMIN) 

nx*-ox 
to  continuf 

ZFRO=l0.n-X«IN)/OX 

print  4 »ox » xm I n ’Zero 
print  4 »0X f xmi n »zfro 

Sl^UxIS  (0..„YSlZr.XRC0.NX.XM2F.0.0.XMIN.DX> 


PRInT  4 ‘OX » XMI n «zfpo 
label  Y axis 


riY=  (ymax-ymini/ysize 
TF  iriY.NF«^'«''M  E’O  To  lA 

FALL  liCALF  (YBFFP«N*Y-SIZF.Yf'lN«nY.ll 
F-O  TO  lA 
U CONTINUF 

no  16  I=1»N 

YbFFR< I )=<XPFFRI I I'YMiNI/nY 

lA  CONTINUF 

print  U.OYtYMIN 

" crt'Ix'.r. 

ny=-ny 

print  lltnYtYMIN 

plot  curvf 

NUMBERI =NUMBFR 

IF  INUMBER*LT»f- ) MUMBrR  = -NH  ‘PE'^ 
Nr)=(N/NUMBFR)-l 

N2  = ° o 

IF  (NUMBER. ME. 1 ) CiO  TO  fl 
IFRST(1)=1 
ILST(1)-N 

« continue 
no  j=i .number 

IF  ( NUMBfRl  .LT.U)  F,0  TT  14fl 
M\  S I E'R&T  ( J ) 

N2=ILST ! j) 

00  TO  1A9 
t4«  CONTINUE 
Nl«NZ+l 
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) 

I 


I 


I 


' * 


K2«Nl*Nn 
1<,0  CONTIfJur 
1PN»1 

tall  plot  IXpFFR(l).YBFFRlnO) 


on  lOO 

print  *i«i5S.XPFFRl  I ) .YPFFRI  J I • 1 
XPF.  XBFFRII) 

YBF.  Yn.'FRt  I ) 

IF  (NSVF.LT.i  ) XPFFRt  n«XSUF-Xp,rFR<  I I 
IF  (XRF.LT.U.O.OPtXpF.f.T.XMZn  f,0  TO  1^0 
IF  ( YPF.lT.u.u.oR.YrF.OT.YSIZF)  GO  TO  lOO 
CALL  PLOT(XnF*YBFtlDM 


lOn 

]«.n 


7S7 


IPN*? 
rONT INUF 
CONTINUF 

IF  I IL.ro. GO  TO  ?5S 

ICT.O 

ICT=1CT*1 

IF  (1CT.GT.1U  GO  To  255  ^ 

CALL  SYMnOL<ILX(  ICT|.1LY(  ICTl  .lLHlICT).lLABLEa»ICT)  .ILRI  ICT) 

1 • ILCI ICT  1 ) 

GO  TO  2R7 


?5A 


?R0 


CONTINUF 
CALL  PLOT 
PRINT  2<)P 

format  (• 


( XSIZF  + R. tO.U.-T) 

plot  noNF*) 


r rfgtorf  the  oata 

nxsARS(nxi 

OY.ARGiOVI 

no  40U  i»i#M 

XPFFRi I l«(OX*XPFFR( I I l+XMIN 
YPFFRI  I ) = (nY»YPFFP(  I ) )*YMIN 

i*oo  continue 

C 

XWAXsO.O 

XMiN'ir.'J 

YMAXsO.O 

Y>'INs0.o 

r 

return 

FNO 

GUGROUTinF  SCALE  ( X ,N « S. YM I N . OY *lO 
OIMFNSIOm  X(2) 

YMAXsXI 1 ) 

ymin=ymax 

NP«N-IC 

00  lU  I »i  *np*k. 

IF( YMAX-X( I ) )5«6*A 
YMAX»X(  I I 

A IF(X(I)-YMIN1  T.lUtlL' 

7 YMIN*X( I I 
10  CONTINUF 

OYa ( YMAX-YHIN) / <S) 
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DO  20  I-l.N 

47300 

?0 

XI  I |.IX( I )-VMIN|/nY 

47400 

print  2S*YMIN«YMAX«dY »X1 1 I »X(mj 

47500 

F0RMATtSE2S.S) 

47600 

return 

47700 

END 

47800 

SUBROUTimE  axis  (X*Y.BCD»NC»SUF»ThETA»YMIN*DY) 

47900 

X1A«,07 

48000 

Slf,N»l  .0 

48100 

IF  (NO  l*?i2 

48200 

1 

sif,N«-i,n 

48300 

? 

NAC«XabSf(NO 

48400 

Th»ThFTA*0,U17A5S2RA 

48500 

NxSl2Ew«S0 

48600 

CTH*C0SF( ThI 

48700 

STH»SIN( TH) 

48800 

TN«N 

48900 

XpxX 

49000 

YB»Y 

49100 

XA  = X-,l*SIfiM*STH 

49200 

YA  = Y+,l»Sir,N*CTH 

49300 

IF  (ThFTa-RU. 111.12*12 

49400 

1? 

XA«-XA 

49500 

in 

CALL  PLOT  (XA.YA.-n 

49600 

49700 

■'f!  LOOP  draws  axis  AMD  MAKES  HC  MARKS 

49800 

49900 

DO  20  I*1.M 

50000 

CALL  PLOT  (XB.YR.21 

50100 

XC»Xn*CTH 

50200 

YCYR+STh 

50300 

CALL  PLOT  (XC.YC.2) 

50400 

• 

XA=XA+CTh 

50500 

ya=ya+sth 

50600 

CALL  PLOT  (XA.YA.21 

50700 

XB  = XC 

50800 

?n 

YB»YC 

50900 

CHARsABSF  (YMINI 

51000 

absv=  absf  (ymin+dY) 

51100 

IF  (ARSV-CHARI  S.6.6 

51200 

ABSVsCHAR 

51300 

A 

FXP=0.0 

51400 

no 

MT=APSV*0.1 

51500 

char=nt 

51600 

IF  (CHAR-ARSVl  <JU. 92*92 

51700 

on 

FXPrEXP-l.r 

51800 

AB5>V«ABiV*10, 

51900 

C,0  To  89 

52000 

'’I 

FXP=FXP+1.0 

52100 

? 

Ap,SV  = ABSV/lO,0 

52200 

NT«ABSV 

52300 

CHAR=nT 

52400 

IF  (CHAR-ApSV19-^.91.91 

52500 

IF  ( ABSV.OT.l  1 C.0  TO  91 

52600 

0-1 

AOY«dY*iO,U**( -EXP  1 

52700 

ARSV=YMn*1L,0**(-fXP1-.TN*ADY 

52800 

XA  = XB-(  .20*SlC.N-0.ns)*STH-.25*CTH 

52900 

YA»YB+(  .2f‘*SlCiN-.0s  j*CTh-.2S*STH 

53000 

N«N41 

53100 

31 


'THET*-90,J  11«3Z«32 
3>  XA»-XA^«1 

no  LOOP  plots  numbers  alomo  axis 
•<1  no  30  !■!  <N 

CALL  NUMBER  1 Xa • Ya • ,07 • AflSVtTHFTA**HF7. 3 » 
ABSV.ABSV-AOY 

IE  UBSV.lt, O.ANO.ABSV.f.T,  11  53.96 
53  ARSV«0 
36  XA»XA-CTh 
30  Ya«YA-STh 
TNC-NAC4.7 

CAtCULATrS  CENTER  POSITION  EOR  LABEL 

XA.X+(SlZr/?,0-,06*TNCl«CTH-(-,C7*S10M*,36)*STH 
YA«Y4.(SlZr/?«0-,06*TNCl*STH-.(-,ti7*SlC.N»,36>»CTM 
IF  (THETa-90,1  41.47./'? 

47  XA>-XA4'.l 

41  CALL  SYMpOL  (XA.VA.X14.nCn.THFT/. MAC) 

PRINTS  LABEL 

XA»XA-.(  (TNr-3,<'l*0,12)«cTH 
YA»Ya.(  (TNC-6.'j|*0,121*STh 
IF  (FXPl  35.30.33 

34  CALL  SYmool  lXA.YA.Xl4.fH(X10  1. THETA. B1 
XAsXA+,4fl*CTH-,07»STH 
YA«YA-.,4B*STH'.,04»CTH 
IF  (EXPl  4".50.4U 

prints  scale  factor 

40  CALI.  NU'.PFR  (XA.YA  .X14.FXP  >ThFTA.4HF3»0| 

30  continue 

CALL  PLOT  (XA.YA. 31 

return 

ENO 

• I CiC. 

' *00 

• • f,r. 
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3-D  Version 


PROGRAM  PHASHIFT 

c Thrff  dimensional  version 

c This  program  calculates  phase  shifts  in  radians  or  degrees  vs  E*A**2 

C OR  K OR  E IN  DIMENSIONLESS  UNITS  FOR  THE  OUANTUM  MECHANICAL  SCATTERING 
C OF  A PArTICLF  of  ENFRGY  E FROM  A ^-0  SOUARt  WELL  OF  WIDTH  A AND  DEPTH 
C “V, 

C The  1ST  data  card  contains  the  following  parameters... 
c The  input  parameters  are  ‘ 1».  e^cro  IThE  lowest  energy  value)  co 

COLS  1-10.(21.  EPSILONIThe  ENERGY  INCREMENT)  COLSll-20.  (?).  SK 1 )»  TO  A 
C (a).  FERMI  (The  well  width)  cols  21-TO.  (5).  MAX  (THE 

C HIGHEST  order  PHASE  SHIFT  DESIRED.  FOR  S.P.DtF)  (1AX=1.2*3»A.  COL  AS 
C (f).  SKIP  To  T (7).  LAST(ThC  TOT 

C AL  number  of  oOINTS  DFSIRED  TO  PE  CALCULATED).  MAX.  1000.  COLS  S5-5B 
C (A)«  JELLO  (IF  printout  of  PHASE  SHIFTS  IN  DFGREES  IS  OESIREOt  ENTER 
C 1 IN  COL  S0«  IF  RADIANS  DESIRFD.  LEAVE  COL  SO  PLANK)  CCLSO 
C ENTER  The  NUmRFR  of  curves  OT  PF  PLOTTFD  IM  COL  60(MAX  9). 

C IF  NO  plot  )S  DESIRED*  ENTFR  A 1 !N  Cf'L.  6S.  IF  ORDINAL  SELF-SCALING 
C DESIRED*  ENTfR  A 1 IN  COL  7t.'«  AND  INSFRT  VALJFS  OF  YMIN  AND  YMAX  REFO 
C T CALL. 

C THF  2ND  data  card  CCNTAINS  THF  WELL  DEPTH  VULCAN*  COLS  1-10.  AND  KILL 
C COL  21 .DESIGNATFS  THE  LAST  DATA  CARD. 

C0MM0N/LARFL/IL*ILX(9) *ILY(o) .IlH(9) .iLRfP) *ILC(9) .ILARLF(9*9) 

common  /range/  XMAX,XMln*YMAX*Vi;N*IFRST(9) *IlST(9)  .xsue 
real  ILX.ILY.ILH* ILR 

real  n 
real  N2 

COMMON/ 1/  delta ( 1000) .CASTOR tlOuu) .POLLUX ( 1 000 ) .HERCULES ( 1 OUO ) .AJA 
IX  ( 1000)  .OLYMPUS!  1000)  .HADES 1 100(i ) .p»  1000  ) .DfU  lOOO)  .))(  1000)  .DN(  1000 
2).B2(lOOO)  .OBZ  ( lOOO  ) »NZ  ( 1 CO.., ) .GAMMA  ( lOOO  ) .ENERGY  ( luOC  ) .CYCLOPS)  100 

30)  .CENT AUR( 1000) .ZFUS( 1000) .RHO I lOOO) .DNZ( 1000) 

DIMENSION  DElTAI ( 1 OoO ) .DELTA? ( 1 OCO ) .DELTA  1000)  .DELTA 

1A( 1000) .W(5) 

dimension  ERGO(lOOO) 

COMMON/3/  PHASEIPOOC) .EROON(ROCO) 

700  read  (60. 800)EZFRO.EPSILON. FERMI .MAX .JELLO  .LAST .NCUHVES. NOPLOT .NOS 

BOO  format ( 3 ( FI  0.7 ) *1 AX » 1 1 »aX . 1 1 *4X . I a . 1 X . 1 1 »aX . 1 1 . AX . 1 1 ) 

0 READ«60.10)VULCAN.KILL 
10  FORMAT(F10,S.10X. 1 1 ) 

IFRSTd)-! 

ifrst(2)«i+  last 

IFRST(3)»1+  2*LAST 
lFRST(A)»l-f  3*LaST 
IFRST(5)«1*  A*LAST 

IFRSTISl-l-t’S'LAST 

IFRST(7)«1*6*LAST 

IFRST(0)»1+7*LAST 

IFRSTI9)»1+8*LA5T 

ilst(i)»  last 

ILST(2)«2*LAST 

ILST(3)«3*LA5T 

ILST(a)»A*LAST 

ILST(5l«5*LAST 

ILST(6)«6*LAST 

TLST(7)"7*LAST 

ILST(8)»8*LAST 

1LST(9)-9*LAST 

XMIN«0. 
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XMAX1»EPSIL0H»LAST  6000 

XTlCK-XMAXl/lO,  6100 

XMAX«<NCURVES-H*XTICX/?.4XMAX1  6200 

XSIZE-XMaX/XTICK  6300 

KC-KC4X  6600 

0RShIFTX-XTICK/2«  6500 

YMIN-0*  6600 

YmAXI-1,  6700 

YTICK-VMAXl/lO,  6600 

ORShIFTY-YTICK/A.  6900 

YMAX««NCURVES-H*0RSHtFTY4YMAXl  7000 

IF  « JELLO.EO.’.J  60  TO  2 7100 

1 PRINT  6 7200 

GO  TO  3 7300 

? print  7 7600 

■>  print  Ip.VULCANtFERMi  »E2ER0»*HAX1  7500 

print  8 7600 

6 format  (iX**PHA6E  ^HIFTS  in  radians  for  a OUANTUh  mechanical  parti  7700 

ICLE  scattering  OFF  A 3-0  SQUARE  WELL»«/)  7800 

7 format  (ix#*phase  Shifts  in  degrees  for  a quantum  mechanical  parti  7900 

ICLE  SCATTFRING  off  a 3-d  SQUARE  WELL.*/)  8000 

R F0RMAT<5Xt*FNERGY( 1 )*6X*0FLTA«0»*6X*DELTa« H*6X*DELTA«2»*6X*DELTA(  8100 
n>*6X*0ELTA«6J*6X*  K *6X*  ERGO  •6X*0LYMPU&1*6X*HADESi*/ ) 8200 

1?  format  I 1X»*WELL  depth-  *F10,5»»  WELL  WIDTH-  *FlO,7**  E^ERO-  * 8300 

1 Fl0,7«»  ELAST-  *F10,7/)  8600 

K-MAX41  8500 

ENERGY!  U-EZFRO  8600 

DO  91  I-1»LAST  8700 

FRGO!  H-FNFRGY!  n»  (FERMI  1**2  8800 

U ENERGY!  UU-EPSILON+EMFRGY!  n 8900 

16  CYCLOPS!  n-SORT!ENERC-Y!  I ) t 9000 

20  CENTAUR!  n«SORT!ENER*'YI  n+VULCAN)  9100 

ZEUS! 1 l-CYCLOPS! n*FERMI  9200 

ZEUSS-ZEUS!  1 1 9300 

27  IF!ZFUSS,F0«0J  78*30  9600 

30  RHO! I 1-CENTAUR! n*FERMI  9500 

RH00*RH0!I)  9S00 

36  CALL  SPHJN!RH00»K*P)  9700 

call  SPHYN(RHOO.K.N)  9800 

DRm»-B(2)  9900 

nNll)«-N!2l  10000 

DO  35  JI-7*MAX  10100 

OBUn-l  JI*BIJI-1»-!JI*1  »*B(  Jl*in/!2*JI*1  ) 10200 

DN<  Jl)-‘ JI*N!  JI-D-!  JI  + l)*N(JI-*in/l2*Jl4l  ) 10300 

7*  continue  10600 

CALL  SPHJN(ZEUSS»K»b2)  10500 

CALL  SPHYN!ZEUSS»K*NZI  10600 

DBZ!U— fli*2)  10700 

DN^Ill— N2!2)  10800 

37  DO  38  KI-2*MAX  10900 

DBi!KII-!KI*RZ(Kl-l  )-(  K I + U*BZ  <KI*1 ))  / ! 2*tCl4U  llOOO 

DMZ!ici|-!Kl*NZ!KI-l)-!ia  + l»*N2(Kl4in/!2*KI  + H lllOO 

3*  CONTINUE  11200 

60  DO  76  LI*1*MAX  11300 

66  GAMMaILII-CFNTaUR!  1 1 *DR(LH/P«Ln  11600 

66  CaSTORIli l-CYCLOPS! n*nBZ!Ll)  11500 

50  POLLUX!lH-GAMMA!LI  I*pZ(Ln  11600 

56  HFRCULFS! LI) -CYCLOPS ! H*0NZ! LI  1 11700 

60  AJAX!LI)-GAHKA!LI)*NZILII  11800 
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OL'^MPUSIl  I ) ■C*STOR<  LI  1-P0LLUX(LI  > 

70  HXPESILI  l.HFRCULESlLn-AJAXtLlt 

7*.  VMLH  -ATAN?tnLYMPUS(Ln  •HAOFSlLin 
IFIWILI)  .LT.O)  VML  I )»WlL  I 
IFI  JFLLF.FO.Ol  r.O  TO  76 
WILI  )»'»tLl  »*A7.2P57fl 
•»(.  roNTINUF 

1F(  JFLLO.FQ.nCiO  To  80 

W(LI)»W<LI  W'.lAlAg 

71  or  TO  bo 

This  sfts  oflta  ■ zero  for  e«o 

7P  00  7P  LI«1 »MAX 
W(L I )«0. 

70  continue 

OLYMPUS! 1 )»0. 

flC  PR^NT  BB .FNFPfiY! I).W(11*WI2I»WIB)»W(6)*W(BI »ZeUS ( I ) »ERGo( ! 1 *0L 
lYMPUS! 1 I .HAOESI 1 1 
n«  format  ciO(?x.Fn.7)  1 

0FLTAln«W(11 
OFLTM  ( iJaW!  ?) 

0FLTA2 ( n *w( 31 
0ELTA3I  n=W(4l 
OFLTAa!  n»wt5) 

01  continue 
BA  print  87  ■ 

SupROUTINE*mIm}maX  automatically  FINOS  YMIN  ANO  YMAX  FOR  LAYOUT. 
Y«IN«0FLTA<MIM ♦ YmaXoOELTa IMAXI f xmin»f2fro*xmax«  EPSIL0N*LAST, 

‘ IF (NOPLOT.FQ.n  00  TO  9S 
IF  (NOSCALF.FO.nOO  TO  99 
CALL  MINIMaX(OELTA«YMIM»YMAX*LAST} 

TMr  plot  routine.  LAYOUT.  IS  DUE  To  STEVE  BRENNER  .CODE  »NRL 
XSiZEiPhYSICAL  length  of  plot  in  inches.  XMIN  =MIN.  VALUE  OF  ABSCISSA 
x*</x»MAx  valuf  of  apscissa.  (xmax-xmin) zxsize  dets  value  of  Tic 
'apJs.  simuIpIy  y refers  to  oroinate.  YSIZE  is  fixed  AT  10  inches. 
In  Th'-'  call  for  subroutine  layout,  the  parameters  in  The  call  state- 
mpkt  are  respectively  apscissa  variable,  ordinate  variable,  number  of 

OaTA  POImTS(=LAST) .no.  of  HOLLERITH  CHARACTERS  FOLLOWED  BY  ABSCISSA 
I^rll.^OrlniNATr  IOFNT.  r.C...  SHOFlTA-O.  (8  characters  MAX).  FOLLOWED 
BY  ThF  actual  number  of  characters  used  in  the  ABSS..ORD.  N*^BS. 
FINALLY  The  Total  nu**rer  of  plots  desired  superimposed  on  Che  Same 

r,pAPH.  IF  CNF  PLOT  PFR  GRAPH.  USE  1.  IP  N.  USE  N.  HOWEVER.  IF  N 

Plots  are  superimposed,  the  the  dimension  statements  for  the  ^rd.abs 
and  also  last  must  be  scaled  by  n.  if  ^clf-scaling  is  desired,  set 

SUBROUTINE  MIMIMAX  AUTOMATICALLY  FINOS  YMIN  AND  YMAX  FOR  LAYOUT. 

Yi<1fi  = DELTA(MlM)  » YMAX  = DEL  TA  I MAX  ) . XM I N*E^ERO  » XMAX«  EPSIL0N*LAST, 

TuF  plot  RCLITinF  WILL  SCALE  THE  PLOTTED  VALUES  PY  POWERS  OF  TEN. 

IF  SriF-SCALlNG  not  USFD8  THEN  YMIN  AND  YMAX  VALUES  MUST  BE  SUPPLIED. 
IF  rrSIRFO.  different  order  PHASF  SHIFTS  MAY  BE  PLOTTED  ON  THE  SAME 
'ooaPhI  SsF  MOSCALF  TO  BYPASS.  PUT  A 1 IN  COL  70.  SUPPLY  YMIN. YMAX. 

no  continue 

NLAST=NCURVFS*LAST 


J«LAST 

ORSHIFTX.npSHlFTX*lcr 

orShifty«orshifty*kc 

jazz«kc*j 

DO  90  M-l.J 

phase (M+JAZZI»DELTA(M)+0RSHIFTY 


11900 

12000 

12100 

12200 

12300 

12400 

12500 

12600 

12700 

12800 

12900 

13000 

13100 

13200 

13300 

13400 

13500 

13600 

13700 

13800 

13900 

14000 

14100 

14200 

14300 

14400 

14500 

14600 

14700 

14800 

14900 

15000 

15100 

15200 

15300 

15400 

15500 

15600 

15700 

15800 

15900 

16000 

16100 

16200 

16300 

16400 

16500 

16600 

16700 

16800 

16900 

17000 

17100 

17200 

17300 

17400 

17500 

17600 

17700 
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FRG0N»M*JAZ2»»ERG0«M)+PRSHIFTX 
RO  CONTInUF 
* PRINT  100 
lOo  format (///j 
RA  IFIKILL.EO.OJ  GO  TO  9 
lFlNOPL.OT,EQ.l)riO  To  199 

CALL  layout ( ERGON *PHASE i NLAST *6HE*A*»2 »7H-0ELTA- #6 • 7 tNCURVES ) 

CALL  STopPLOT 
IRp  END 

subroutine  minimaxi  array tYMlN»VMAX. LAST) 
dimension  ARRAYIlOOO) 

c This  subrout  inf  finds  both  the  largest  and  smallest  numbers  in  array, 

C LARGEST  NUMBER 
■JJ*1+J 

10  DO  30  K.JJtLAST 

TF I ARRAYIJ).GF.ARRAYIK 1)30*20 
20  J«X 

IFlK.EQ.LAST)  GO  TO  40 
GO  TO  30 
30  CONTINUF 
40  YMAX*ARRAY(J) 

PRINT  43. JJ 
4A  FORMAT  (1X.*JJ«*I3) 

GO  TO  60 
30  JJ»J+1 


GO  irO  10  20400 

C smallest  NUMRFR  20500 

60  L»1  20600 

LL=1+L  20700 

lOO  DO  300  M«LL*LAST  20800 

IF<ARRAY(L).LF«ARRAY(M) ) 300*200  20900 

2C0  L»M  21000 

IF<M,EQ.UAST)G0  to  400  21100 

GO  TO  500  21200 

3OO  CONTINUE  21300 

4O0  Ymin=ARRAY(L)  21400 

print  430*LL  21500 

430  FORMAT  IlX**LL=«I3)  21600 

GO  TO  600  21700 

3O0  LL=L+1  21800 

GO  TO  lOO  21900 

6O0  return  22000 

END  22100 

subroutine  SPHJNIX*n*ARR)  . 22200 

C3  NpL  RESN  BESSEL  FUNCTIONS  OF  THE  FIRST  AND  SECOND  KINDS  22300 

DIMENSION  ARRIll  22400 

RXsi,0E+307  22500 

NF*1  22600 

NN»N+1  22700 

IFtX.EO.O.C)  GO  TO  4 22800 

IFIX,LE.0,5)G0  to  6 22900 

C FOR  ORDERS  LESS  THAN  ARGUMENT*  AND  ARGUMENT  GREATER  THAN  0,5  23000 

APR(1)>SiN(X)/X  23100 

ARR(2)=ARR(1)/X-C0SIX)/X  23200 

NMaX-2  23300 

IF(N,lT,X)NM«N-1  23400 

IFINM.LT,0)NM»0  23500 

DO  16  I-l.NM  23600 


17800 

17900 

18000 

18100 

18200 

18300 

18400 

18500 

18600 

18700 

18800 

18900 

19000 

19100 

19200 

19300 

19400 

19500 

19600 

19700 

19800 

19900 

20000 

20100 

20200 


t 


FN-I+I  23700 

ARR<  I-f2).(FN4l.0)»ARR<  I + l )/X-ARR(  n 23800 

C CHECK  FOR  UNDERFLOW  23900 

IF«ARRI  I'f2i«NE.0.0)GO  TO  16  24000 

ARRI  l4’l)>ARR(I'fl)»RX  24100 

ARR(  l4-2)«lFN'fl.U)*ARR(  I-fl  l/X-ARR(  I)«RX  24200 

PRINT  13(  XdtN  24300 

16  continue  24400 

IF«N.LT,X)RFTURN  24500 

H0L0*ARRINM'«-2I  24600 

C FOR  ORDERS  GREATER  THAN  OR  EQUAL  TO  ARGUMENT*  AND  ARGUMENT  24700 

C greater  than  0.5  24800 

lA  ALPhA*1.0  24900 


M.-2-NM 

NOTE-O 

I»N+47 

IFfX/N.GT.O.on-I+X/li 

ARRI l42)aALPHA 

ARRI 1+31. 0,0 

I«-I-l 

DO  2 XK»I*M 

X»-KX 

N00»K-*1 

TxPl.X-fX  + 1 

ARR(KI*T)cP1*ARR(K-*1  )/ 
IF  exponent  fault  18* 
1«  TeMPY.ARRIic  )«ARR<K+1 '( 
ARR(X)»ARR(K)*RX 
IF  EXPONENT  FAULT  24* 
2fc  ARR(X>»TeMPY 


ARRlX+llaARRlX+l) /RX  26700 

• print  13*X*M*N00  26800 

? continue  26900 

12  C«HOLD/ARR(r.)  27000 

IF(C.NE.O.O)GO  TO  U 27100 

, H0L0»H0LD*RX  27200 

N0TE»N0Tr*-l  27300 

GO  TO  1?  27400 

U IFINOTE.EO.OIGO  TO  17  27500 

HULT=107*NOTE  ' 276C0 

NUL=L-1  27700 

PRINT  21 *X*N*NUL*MUlT  27800 

71  format  (3X*F0R  ARGU'-IENT  OF  *E17»lO»  AND  ORDERS  OF  *14*  TO  *14*  27900 

1 *SPHJN  exponent  was  UJCREASED  BY  *16)  28000 

17  DO  3 I=IC*NN  28100 

ARR(  II*ARR(  n*C  28200 

3 continue  28300 

I return  28400 

I A ARR«1)=1,0  28500 

i DO  5 L»2*NN  28600 

I ARR«L)»0,0  28700 

I A continue  28800 

return  28900 

c series  for  argument  less  Than  or  eoual  to  0,5  29000 

A FPS.5.OF-II  29100 

CH*CS«DN«P1V»1,0  29200 

DO  lO  L«1*NN  29300 

. ITFXP*NEXPA«MEXPA«LFXPA"0  29400 

NC»L-1  29500 


1 


X-ARRIK+2) 

2 

/RX*TkP1/X-ARRIK*2)  /RX 


25000 

25100 

25200 

25300 

25400 

25500 

25600 

25700 

25800 

25900 

26000 

26100 

►26200 

26300 

26400 

26500 

26600 
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IFtL'.EO.l)  GO  TO  7 
CS»CS*X/pN 
CHECK  FOR  underflow 
IFCCS.NE.0,0)  GO  TO  7 

Scale  for  this  n and  all  following  n 

CS»CH*RX*X/DN 
NEXPA-NEXPA+1 
7 DN-DN+2.0 
CH-CS 
XKK-0.0 

A»CH2»FOrTT«1,0/OIV 
Tn»2.0*NO 
FN«-0.5*X*X 
n xkk-xkk+1.0 

FD»XKK*(TN+XKK+XkK>1.0» 

A.A*FN/FD  .*  • 

CHECK  FOR  underflow  • . 

IFU.NE.O.OIGO  To  25 

Scale  for  this  n and  a\l, following  n 

A»CH2*RX*FN/F0  ' . ’ 

div.div/rx  • 

MFXPAiMEXPAfl  • 

FORTT-FORTT*RX  ’ ' 

25  CH?«k  • 

TTmA+FORTT  , * 

IF-<TJ,,EQ.O.O,OR,APS(A/TTl,LT.EPSl  GO  To  9 

FORTTsTT 

GO* TO  8 

o ARR(L)«CS*TT 

CHFCK  FOP  UNDERFLOW 

IF( ARR(LI .NF.O.OIGO  To  31 

scale  for  THIS  N AND  ALL  FOLLOWING  N 

ARR(H»CS*RX*TT 

DIV»DIV/RX 

LEXPA«LEXPA*1 

31  itexp=nfxpa*mexpa+lexpa 
IF  I lTEXP.ro. CJGO  TO  10 
IF( ITEXP.F0.1IG0  TO  29 
PRINT  30.  NEXPAtMEXPAfLEXPA 
30  F0RmaT(3X.»  SCALF  FRROR  *3<3X»U0)I 
?0  print  13*  X.N0*N 

!■>  F0RMAT«3X.*F0R  argument  of  **E17.10»*  and  orders  of  **U*«  To  *• 
1 14.*  SPHJf!  EXPONENT  WAS  INCREASED  BY  307*) 

10  continue 

END 

subroutine  SphYnI X.N.ARR) 
dimension  ARR(1) 

RX.I.Oe+307 
IF<X.EO.O.)5*3 
A print  4 

4 format <3X. *00  NOT  USE  ZERO  ARGUMENT  FOR  SPHYN*) 

return 

FIND  YO  AND  Y1 
3 ARR( 1 )«-COSIX)/X 

ARR(2)>ARR( 1)/X-SIN<X)/X 
NMlmN-1 
DO  2 K*1»NM1 
TKP1«K+K4.1 

ARR(K*2)»TKP1*ARRIK*1)/X-ARR«K) 


29600 

29700 

29800 

29900 

30000 

30100 

30200 

30300 

30400 

30500 

30600 

30700 

30800 

30900 

31000 

31100 

31200 

31300 

31400 

31500 

31600 

31700 

31800 

.31900 

32000 

32100 

32200 

32300 

32400 

32500 

32600 

32700 

32800 

32900 

33000 

33100 

33200 

33300 

33400 

335C0 

33600 

33700 

33800 

33900 

34000 

34100 

34200 

34300 

34400 

34500 

34600 

34700 

34800 

34900 

35000 

35100 

35200 

35300 

35400 
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C CHECK  FOR  OVERFLOW  35500 

IF  exponent  fault  8t2  35600 

« ARRdCvlIaARRIK't'D/RX  35700- 

ARR<K-f2)«TKPi*ARR(K-»n/X-ARRUj/RX  35800 

print  7»X*KtN  35900 

7 F0RMAT(3X«*F0R  argument  of  *»F17«1C»*  and  orders  of  *»I*t*  To  *»  36000 

1 1<,,*  SPhVN  exponent  was  decreased  by  307*)  36100 

2 continue  36200 

END  36300 

subroutine  layout  (XBFFR.YBFFR.NN#XBCD*YBC0*NNX»NNY. number)  36400 

COKMON/LAOEL/ILt ILX(9) 1 1 L Y«9 1 ♦ ILH« 9 ) • ILR 1 9 1 » I LC • 9 ) * I LABLE « 9*9 ) 36500 

COMMON  /RANGE/  XvaX tXM IN t YMAX » YMIN» IFRST ( 9 ) * ILST ( 9 ) iXSIZE  36600 

real  ILX.ILYtILH»ILR  3670(f 

dimension  XBCDI5) .YpCD<51  36800 

dimension  XRFFPiNNt «YRFFR(NN> »PLTAPRAY<254»  36900 

DATa  (lSTRT«n  37000 

C initialization  37100 

H»NN  37200 

NX=nnX  37300 

NY=nNY  37400 

IF  (XSIZf.FO.0,0)  XSIZE-IO.  37500 

YSI2C«1<-.  37600 

NSVE=N  37700 

N«XabSF  (N)  37800 

IF  (number. FO.OI  NUMBFR»1  37900 

IF  (ISTRT.FO.2)  GO  TO  2 38000 

CALL  PLOTS(PlTARRAY.254.131  38100 

ISTpT=2  38200 

? CCNTINUF  38300 

38400 

LAPFL  X AXIS  38500 

38600 
38700 

nXx  (XMAX-XmINI/XSIZE  38800 

IF  (PX.NF.U)  GO  TO  3 38900 

CALL  scale  (XRFFR.N»XSIZE*XMIN»DX.1 ) 39000 

GO  TO  7 39100 

t CCNTINUF  39200 

DO  6 I«1.N  39300 

r PRINT  5S5S.XPFFR(  n .YPFFR(n»I_  39400 

«i5'!e  F0R'»AT  I?F70,lo.IlO)  ' 39500 

A XBFFR(I)=  (XRFFRI n-XMlN)/nX  39600 

7 CONTINUF  39700 

39800 

39900 

C CALL  SCAlF  (XBFFR.NfXSlZE.XMIN.DX.l ) 40000 

ZeR0=XSIZF-ZER0  40100 

IF  (NSVEJ  5*10.10  40200 

A Xvin«XSIZF*DX-ABSF( XMlNl  4O3OO 

DX«-DX  40400 

in  CCNTINUF  40500 

2ER0»(0,0-XfIN)/DX  40600 

print  A*DX*XMIN*ZFR0  40700 

4 FCRMAT  I*  DX«**3F20,51  40800 

CALL  AXIS  (0,0.0, O.XBCD*NX»XSIZF»v',O.XMIN.DX)  40900 

print  4«DX,XMIfi.ZER0  41000 

NX. -NX  41100 

CALL  AXIS  (i;,0,YSlZF*XBCO.NX.XSIZE.O,O.XMlN.DX)  41200 

r 41300 
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PRINT  4ipX«XMINtZrP0 
UAREL  T AXIS 


t)Y.  IVMAX-VMINI/VSI2E 
IF  (OY.NE.U.OI  go  to  14 
CALL  bCALC  I YflFFR.N . YS I ZE • YMI N*OV .1 ) 
GO  TO  15 

14  continue 

00  16  1»1*N 

U YBFFRm-IYRFFR(I)-YMlNl/r>Y 

15  CONTINUE 


‘ print  11.0''.YMIN 
ll  format  (•  0Y»*«3F2^»4» 

CALL  AXIS  IO,O.O.O.YBCD»NY«YSlZFt90.»YMlN»DY) 
NY. -NY 

PRINT  lltfjY.YMlN 

CALL  AXIS  (XSlZE*0.f'.YBC0«NY,YSIZE*90..YMIN.DY) 
PRINT  lltDY.YMIN 

PLOT  CURVE 
NUMBFR1»NUMBFR 

IF  t number. LT, 01  NUMBER«-NUMBER 

NC«(N/NUMPERJ-1 

N?»0 

IF  (NUMPFR.ME.1 ) GO  TO  B 

IFRSTd  1.1 

ILSTIll.N 

A continue 

00  15U  J.liNUMBER 

IE  tNUMBFRl.LT.O)  GO  TO  148 

N1»IFRST(J) 

N2=1LSTI j) 

GO  To  149 

148  continue 

N1«N2>1 

N2.N1+N0 

140  continue 

IPN.3 

CALL  PLOT  (XBFFRd)  .YBEER(H»3J 
00  lOO  I.N1.N2 

PRINT  5555.XBFFR( I)  .YBFFRt n .1 
XBF«  XBEFRm 
YBF.  YBFFRt I I 

IE  (NSVE.lt. 0)  XBFERI  n»XSUF-XBFFR(  I > 

IE  (XPE.lT.O.O.OR.XbF.GT.XSIZE)  go  to  100 
IF  ( YBF.lT.O.O.OR.YpF.GT.YSIZF)  go  to  lOU 
CALL  PL0T(XBF»YRF.IPN) 

IPN.2 

ICO  CONTINUE 

150  continue 

IF  IIL.FO.0I  GO  TO  255 
ICT-0 

?5?  ICT.ICT*! 


41400 

41500 

41600 

41700 

41600 

41900 

42000 

42100 

42200 

42300 

42400 

42500 

42600 

42700 

42B00 

42900 

43000 

43100 

43200 

43300 

43400 

43500 

43600 

43700 

43800 

43900 

44000 

44100 

44200 

44300 

44400 

44500 

44600 

44700 

44800 

44900 

45000 

45100 

45200 

45300 

45400 

45500 

45600 

45700 

45800 

45900 

46000 

46100 

46200 

46300 

46400 

46500 

46600 

46700 

46800 

46900 

47000 

47100 

47200 


JCT?.  ICT.  MLH(  ICT)  .ILA^^^ 

iMLcncTn 
GO  TO  252 
255  CONTINUF 

CALL  PLOT  <XSlZF+3.*0.0.-3» 

PRINT  250 

250  FORMAT  t*  PLOT  OONE*I 


rfstorf  the  data 

pX.ABSfOX) 

OY.ABSInYt 

00  aOO  I«ltN 

XPFFRi 1 )b(OX*XHFFRI I ) 1+XMlN 
YBFFR<  t I = (OY*YbFFR< I J )+YM!N 
AOO  CONTINUF 

xmax.o.u 

XMIN«0.0 

YMAX=0.0 

YM1N=C.0 

return 

SUPROUTINF  scale  ( X.N»S*YMtNtOY.K) 

OIMENSION  X(7) 

YyAXeXn  ) 

YMIN-YMAX 

NP«N*X 

00  lO  I«1»NP»X 
IF(YMAX-X(I))5*6*5 
S YMAXuXd) 

A IF«Xm-YMIN)  TtlO.lO 
7 YMIN«X( I ) 

in  continue 

rY«( YMAX-YMIN) / (SI 
00  20  I«1*N 

?0  X(  nslX(  I )-YMIN)/0Y 

print  25tYMIMtYMAX»oY*Xll) »XtN) 

25  F0RMAT(5F25.5J 

rfturn 

subroutine  axis  (XtYtBCO»NC»SIZF*THETA*YMlN»DY) 

X1A«,07 

SlGN»l .0 

IF  (NO  1»2»2 

1 SIGN«-1»0 

7 nac«xabsf<no 

Th=THETA*0,017A5329A 
N«StZE+«50 
CTH«COSF( TH) 

STH»SIN(THJ 

Tn«N 

Xfl-X 

YB«Y 

XA«X-,1*SIGN*STh 


41 


r\  r>  ry 


yA»Y+,l*S!riN«CTH 
r IF  ITHETA-90,)lU,12tl2 

C 1>  XA»-XA 

10  CALL  PLOT  (Xa»YA*3) 


DO  LOOP  DRAWS  AXIS  AND  MAXES  HC  MARKS 


DO  20  I»13M 

CALL  PLOT  (Xp.YR,2) 

XC«XBt-CTH 


yc«yb+sth 

CALL  PLOT  (XC.YC.2I 


XA«XA+CTh 

Ya-YA+STh 

CALL  PLOT  (XA.YA.2) 


xp«xc 

20  YRcYC 

CHAR»ABSF  (YMINI 
ARSV«  APSF  (YMIN^DY) 

IF  (ABSV-CHAR)  5*6*6 
A APSV*CHAR 
f FXPxO.O 
RR  NT  = ABSV4.0,1 
CHAR-iNT 

IF  (CHAR-APSV)  90*92*92 
oo  FXP«EXP-1,0 
ARSV=ABSV*10« 

GO  TO  89 
FXP  = FXP-*1,0 
9?  ABJ>V»ABSV/10,0 
NTsAPSV 
char*nt 

IF  (CHAR-ARSV)93*91 *91 
IF  (AHSV.GT.ll  GO  TO  91 
9*  ADYst)Y*lO,ii»*(-FXPl 

ARSVaYMlN»lO,0**( -EXP)+TN*A0Y 
XA=XR-( ,20#S1GN-0,05)*STH-.25*CTH 
Ya=Yp+( ,20*SIGN-,0r j*cTh-,25*STh 
N = N+1 

IF  (ThfTa-90.)  31*32*32 
3?  XA«-XA*,1 


r on  LOOP  PIOTS  mumrfps  along  axis 

c 

31  DO  30  I«1*M 

CALL  number  tXA*YA*,07*ABSV*THETA*AHF7,3) 
ARSVxABSV-ADY 

C IF  (ABSV.LT»o,anD.APSV.GT,-H  55*56 

r A*  Ap,sv=o 
SA  XA=XA-CTh 
30  YAxYA-STh 
TNCxNAC3-7 


c calculates  center  position  for  label 

XA»X+(SiZF/2.0-.06*TNC1*CTH-(-»07+SIGN*,36)*STH 
YAxY+(SlZF/2.0-.06*TriC)  ■•‘oTH+I-,L7+SIGN».36)«CTH 
r IF  (ThETa-90.)  41*A2*A2 

c 42  XA=-XA+,1 

41  call  symbol  IXA*Ya*X14*PCD*THETA*NAO 


' -"1 


t • 


C PRINTS  LABEL 

XA-XA4I  (TNC-f>.'Jl*0.12»*CTH 
Ya»Ya+) ITNC-6.01*0,12J*STh 
IF  (EXP)  35»50.15 

3A  CALL  SYMBOL  ( XA . YA »XI A •OH ( XlO  )»THETAt8) 

XA»XA+.48«CTH-.07*STH 
Ya=YA+.A8*STh+.04*CTH 
IF  (EXP)  4C*50t40 
C 

c prints  scalf  factor 

40  CALL  number  ( XA . Ya  • X 14 »FXP»TheTA .4HF3 .0 ) 

50  CONTINUF 

CALL  PLOT  (XAtYA.3) 

return 

END 

• 'GCi 

• *00 

• • GCi 


59100 

59200 

59300 

59400 

59500 

59600 

59700 

59800 

59900 

60000 

60100 

60200 

60300 

60400 

60500 

60600 


I » 


43 


t 

i • 


